
# THIS FILE CONTAINS CODE FOR RECREATING THE FOLLOWING MATERIALS FOR "Do 
# Politicians Outside the United States Also Think Voters are More Conservative 
# than they Really Are? A Comparative Study of the Conservative Bias in Elites’ 
# Perception of Public Opinion"
#
# THE FOLLOWING ARE CREATED USING THE CODE HERE:
# MANUSCRIPT FIGURE 1 + ACCOMPANYING SI TABLE A4.1
# MANUSCRIPT FIGURE 2 + ACCOMPANYING SI TABLE A4.2
# MANUSCRIPT FIGURE 3
# MANUSCRIPT FIGURE 4
# SI FIGURE A4.1 + ACCOMPANYING SI TABLE A4.3
# SI FIGURE A5.1 + ACCOMPANYING SI TABLE A4.4
# TABLES A1.1, A1.2, A1.3, A2.1, A3.1, A7.1, A8.1 IN SI WERE PRODUCED MANUALLY BY THE AUTHORS. 

library(stargazer)
library(foreign)
library(reshape2)
library(dplyr)
library(tidyr)
library(ggplot2)
library(scales)
library(grid)
library(readstata13)
library(tibble)
library(modelr)
library(psych)
library(RColorBrewer)
library(lme4)
library(lmerTest)

# SETUP
rm(list = ls())
theme_set(theme_grey())
# Fill in with relevant working directory
setwd("[YOUR WORKING DIRECTORY]")

df <- read.csv("pilet_et_al_replication_data_politicians.csv")

View(df)

# DATA PREP
# =========

# Create country names variable
df$country_name <- NA
df$country_name[df$V2==1]<-"Flanders"
df$country_name[df$V2==2]<-"Wallonia"
df$country_name[df$V2==3]<-"Switzerland"
df$country_name[df$V2==4]<-"The Netherlands"
df$country_name[df$V2==5]<-"Germany"
df$country_name[df$V2==6]<-"Canada"

# FILTER OUT NL
df <- df %>%
  filter(country_name!="The Netherlands")  

# party classification variables
df$party_lr_bin <- NA
df$party_lr_bin <- ifelse(df$party_family<=30,"Left",
                          ifelse(df$party_family>30&df$party_family<=95,"Right",NA))
df$party_lr_bin[df$V4=="SVP"]<-"Right"
df$party_lr_bin[df$party_family==999]<-NA

df$party_lcr<-NA
df$party_lcr <- ifelse(df$party_family<=30,"Left",
                       ifelse(df$party_family>50&df$party_family<=95,"Right",
                              ifelse(df$party_family>=40&df$party_family<=50,"Centre",NA)))
df$party_lcr[df$V4=="SVP"]<-"Right"
df$party_lcr[df$party_family==999]<-NA

# leave only statements with clear liberal/conservative direction
df <- df %>%
  filter(conservative_direction!=99)

# politicians' estimations of country public opinion
df$estimation_country_agree<-df$Q12
df$estimation_country_agree[df$estimation_country_agree==-99]<-NA

# politicians' estimations of electorate public opinion
df$estimation_electorate_agree<-df$Q16
df$estimation_electorate_agree[df$estimation_electorate_agree==-99]<-NA

# politicians' esitmations of district PO (Switzerland)
df$estimation_district_agree<-df$Q12_dis
df$estimation_district_agree[df$estimation_district_agree==-99]<-NA

# remove -99 in cou_agr_cal, ele_agr_cal, dis_agr_cal
df$cou_agr_cal[df$cou_agr_cal<0]<-NA
df$ele_agr_cal[df$ele_agr_cal<0]<-NA
df$dis_agr_cal[df$dis_agr_cal<0]<-NA

# create country-level CB variable
df$cb_country_level<-ifelse(df$conservative_direction==1,
                            df$estimation_country_agree-df$cou_agr_cal,NA)
df$cb_country_level<-ifelse(df$conservative_direction==0,
                            (-1)*(df$estimation_country_agree-df$cou_agr_cal),df$cb_country_level)

# create electorate-level CB variable
df$cb_electorate_level<-ifelse(df$conservative_direction==1,
                               df$estimation_electorate_agree-df$ele_agr_cal,NA)
df$cb_electorate_level<-ifelse(df$conservative_direction==0,
                               (-1)*(df$estimation_electorate_agree-df$ele_agr_cal),df$cb_electorate_level)

# create district-level CB variable
df$cb_district_level<-NA
df$cb_district_level<-ifelse(df$conservative_direction==1,
                             df$estimation_district_agree-df$dis_agr_cal,NA)
df$cb_district_level<-ifelse(df$conservative_direction==0,
                             (-1)*(df$estimation_district_agree-df$dis_agr_cal),df$cb_district_level)

# create country and district-level CB variable
df$cb_country_and_district_level<-df$cb_country_level
df$cb_country_and_district_level<-ifelse(!is.na(df$cb_district_level),df$cb_district_level,df$cb_country_level)

# Create by-stataement topic variable
df$topic_name <- NA
df$topic_name[df$SID %in% c("be_5","ch_A8","ch_B3","de_A3","de_B3", "ca_2")] <- "State Intervention"
df$topic_name[df$SID %in% c("be_8","ch_A9","de_A6","ca_6","nl_2")] <- "Pension"
df$topic_name[df$SID %in% c("be_6","ch_A6","de_B8","de_A8","ca_1","nl_3")] <- "Immigration"
df$topic_name[df$SID %in% c("ch_A4","ch_B4", "ca_8")] <- "Cultural"

# 95% ci functions
lower_ci <- function(mean, se, n, conf_level = 0.95){
  lower_ci <- mean - qt(1 - ((1 - conf_level) / 2), n - 1) * se
}
upper_ci <- function(mean, se, n, conf_level = 0.95){
  upper_ci <- mean + qt(1 - ((1 - conf_level) / 2), n - 1) * se
}


# FIGURE 1
# ========

df_cbplot_country_and_district <- df %>%
  group_by(country_name,topic_name) %>%
  summarise(cb_n = sum(!is.na(cb_country_and_district_level)),
            cb_mean = mean(cb_country_and_district_level,na.rm=TRUE),
            cb_ssd = sd(cb_country_and_district_level,na.rm=TRUE)) %>%
  mutate(cb_se = cb_ssd / sqrt(cb_n),
         cb_cilo = lower_ci(cb_mean, cb_se, cb_n),
         cb_cihi = upper_ci(cb_mean, cb_se, cb_n))

df_cbplot_electorate <- df %>%
  group_by(country_name,topic_name) %>%
  summarise(cb_n = sum(!is.na(cb_electorate_level)),
            cb_mean = mean(cb_electorate_level,na.rm=TRUE),
            cb_ssd = sd(cb_electorate_level,na.rm=TRUE)) %>%
  mutate(cb_se = cb_ssd / sqrt(cb_n),
         cb_cilo = lower_ci(cb_mean, cb_se, cb_n),
         cb_cihi = upper_ci(cb_mean, cb_se, cb_n))

df_cbplot_country_and_district$estimation_type="Country / District"
df_cbplot_electorate$estimation_type="Electorate"
df_cbplot_country_and_district_and_electorate<-df_cbplot_country_and_district %>%
  bind_rows(df_cbplot_electorate)

df_cbplot_country_and_district_and_electorate$topic_name[df_cbplot_country_and_district_and_electorate$topic_name=="State Intervention"]<-"Economic"

myColorScale <- brewer.pal(4,"Dark2")
names(myColorScale) <- levels(df_cbplot_country_and_district_and_electorate$topic_name)
cb_countries_by_topic_name_and_level_including_district_facets<-ggplot(data=df_cbplot_country_and_district_and_electorate)+
  geom_pointrange(aes(x=reorder(country_name,desc(country_name)),
                      y=cb_mean,ymin=cb_cilo,ymax=cb_cihi,group=cb_mean,
                      #colour=topic_name,
                      shape=estimation_type),
                  position=position_dodge(width=0.5))+
  geom_hline(aes(yintercept=0),linetype="dashed")+
  facet_wrap(~topic_name)+
  ylab("Conservative Bias, Percentage Points")+
  xlab("")+
  scale_shape_manual(name = "Estimation Type",values=c(16,17))+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=14),
        legend.title = element_text(size=14),
        legend.text = element_text(size=12),
        strip.text.x = element_text(size = 12))+
  coord_flip()
cb_countries_by_topic_name_and_level_including_district_facets
ggsave("cb_countries_by_topic_name_and_level_including_district_facets.pdf",cb_countries_by_topic_name_and_level_including_district_facets)
ggsave("cb_countries_by_topic_name_and_level_including_district_facets.png",cb_countries_by_topic_name_and_level_including_district_facets)

# DATA TABLE FOR FIGURE 1  [SI TABLE A4.1]
# ========================================

df_cbplot_country_and_district_and_electorate %>%
  select(-c(cb_ssd,cb_se)) %>%
  as.data.frame() %>%
  stargazer(summary=FALSE,type="html")

# FIGURE 2
# ========

df_cbplot_issue_party_lr_country_and_district <- df %>%
  group_by(country_name,topic_name,party_lr_bin) %>%
  summarise(cb_n = sum(!is.na(cb_country_and_district_level)),
            cb_mean = mean(cb_country_and_district_level,na.rm=TRUE),
            cb_ssd = sd(cb_country_and_district_level,na.rm=TRUE)) %>%
  mutate(cb_se = cb_ssd / sqrt(cb_n),
         cb_cilo = lower_ci(cb_mean, cb_se, cb_n),
         cb_cihi = upper_ci(cb_mean, cb_se, cb_n))

myColorScale <- brewer.pal(4,"Dark2")
names(myColorScale) <- levels(df_cbplot_issue_party_lr_country_and_district$topic_name)
cb_countries_country_and_district_by_issue_party_lr_facets<-ggplot(data=subset(df_cbplot_issue_party_lr_country_and_district,
                                                                               !is.na(df_cbplot_issue_party_lr_country_and_district$party_lr_bin)))+
  geom_pointrange(aes(x=reorder(country_name,desc(country_name)),
                      y=cb_mean,ymin=cb_cilo,ymax=cb_cihi,group=cb_mean,
                      shape=party_lr_bin),
                  position=position_dodge(width=1))+
  geom_hline(aes(yintercept=0),linetype="dashed")+
  facet_wrap(~topic_name)+
  ylab("Conservative Bias, Percentage Points")+
  xlab("")+
  scale_shape_manual(name = "Party Ideology",values=c(16,17))+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=14),
        legend.title = element_text(size=14),
        legend.text = element_text(size=12),
        strip.text.x = element_text(size = 12))+
  scale_shape_manual(values=c(22,23))+
  guides(shape=guide_legend(title="Party Ideology"))+
  coord_flip()
cb_countries_country_and_district_by_issue_party_lr_facets
ggsave("cb_countries_country_and_district_by_issue_party_lr_facets.pdf",cb_countries_country_and_district_by_issue_party_lr_facets)
ggsave("cb_countries_country_and_district_by_issue_party_lr_facets.png",cb_countries_country_and_district_by_issue_party_lr_facets)


# DATA TABLE FOR FIGURE 2  [SI TABLE A4.2]
# ========================================

subset(df_cbplot_issue_party_lr_country_and_district,
       !is.na(df_cbplot_issue_party_lr_country_and_district$party_lr_bin)) %>%
  select(-c(cb_ssd,cb_se)) %>%
  as.data.frame() %>%
  stargazer(summary=FALSE,type="html")


# FIGURE 3
# ========


# create mean CB country and district level score for each politician
df_politician_country_and_district_mean <- df %>%
  group_by(X) %>%
  summarise(cb_country_and_district_level_politician_mean = mean(cb_country_and_district_level,na.rm=TRUE),
            country_name=min(country_name)) %>%
  ungroup() %>%
  filter(!is.na(cb_country_and_district_level_politician_mean))

View(df_politician_country_and_district_mean)
ggplot(data=df_politician_country_and_district_mean)+geom_histogram(aes(x=cb_country_and_district_level_politician_mean),color="gray53",fill="white",bins=50)+
  geom_vline(aes(xintercept=0),size=0.5,linetype="dashed")+
  facet_wrap(~country_name,scales = "free_y")+
  xlab("")+
  ylab("Number of Politicians")+
  theme(strip.text = element_text(size=12),
        axis.text = element_text(size=12),
        axis.title = element_text(size=12))
ggsave("cb_country_and_district_full_histogram_by_country.pdf")
ggsave("cb_country_and_district_full_histogram_by_country.png")


# FIGURE A5.1 IN SI
# =================

df_cbplot_country_district_issue_party_lcr <- df %>%
  group_by(country_name,topic_name,party_lcr) %>%
  summarise(cb_n = sum(!is.na(cb_country_and_district_level)),
            cb_mean = mean(cb_country_and_district_level,na.rm=TRUE),
            cb_ssd = sd(cb_country_and_district_level,na.rm=TRUE)) %>%
  mutate(cb_se = cb_ssd / sqrt(cb_n),
         cb_cilo = lower_ci(cb_mean, cb_se, cb_n),
         cb_cihi = upper_ci(cb_mean, cb_se, cb_n))


myColorScale <- brewer.pal(4,"Dark2")
names(myColorScale) <- levels(df_cbplot_country_district_issue_party_lcr$topic_name)
cb_countries_country_district_by_issue_party_lcr<-ggplot(data=subset(df_cbplot_country_district_issue_party_lcr,!is.na(df_cbplot_country_district_issue_party_lcr$party_lcr)))+
  geom_pointrange(aes(x=reorder(country_name,desc(country_name)),
                      y=cb_mean,ymin=cb_cilo,ymax=cb_cihi,group=cb_mean,
                      colour=party_lcr,shape=party_lcr),
                  position=position_dodge(width=0.5))+
  geom_hline(aes(yintercept=0),linetype="dashed")+
  facet_wrap(~topic_name)+
  ylab("Conservative Bias, Percentage Points")+
  xlab("")+
  scale_colour_manual(name = "Party Ideology",values = myColorScale)+
  guides(shape="none")+
  theme(axis.text=element_text(size=12),
        axis.title=element_text(size=14),
        legend.title = element_text(size=14),
        legend.text = element_text(size=12),
        strip.text = element_text(size=12))+
  coord_flip()
cb_countries_country_district_by_issue_party_lcr
ggsave("cb_countries_country_district_by_issue_party_lcr.pdf",cb_countries_country_district_by_issue_party_lcr)
ggsave("cb_countries_country_district_by_issue_party_lcr.png",cb_countries_country_district_by_issue_party_lcr)


# DATA TABLE FOR FIGURE A5.1 [SI TABLE A4.3]
# ==========================================

subset(df_cbplot_country_district_issue_party_lcr,
       !is.na(df_cbplot_country_district_issue_party_lcr$party_lcr)) %>%
  select(-c(cb_ssd,cb_se)) %>%
  as.data.frame() %>%
  stargazer(summary=FALSE,type="html")


# Plotting citizen activism by issue and liberal/conservative position
# =====================================================================

df_citizen_activism <- read.csv("polpop_citizen_activism_by_conservatism.csv")

# FIGURE 4
# ========

myColorScale <- brewer.pal(5,"Dark2")
names(myColorScale) <- levels(df_citizen_activism$issue_type)
citizen_activism_politician_contact_plot <- ggplot(data=subset(df_citizen_activism,activity=="Contacted Politicians"))+
  geom_pointrange(aes(x=issue_type,
                      y=mean,ymin=cilo,ymax=cihi,group=mean,
                      colour=issue_type,shape=position_on_issue),
                  position=position_dodge(width=0.5),size=0.8)+
  scale_colour_manual(name = "Policy Category",values = myColorScale)+
  scale_shape_manual(name = "Position on Issue",values = c(16,17))+
  xlab("")+
  ylab("Mean Share of Citzens Who Contacted a Politician")+
  theme(axis.text.x = element_text(size=12,angle=30,hjust=1),
        axis.text.y = element_text(size=12),
        axis.title.y = element_text(size=12),
        legend.text = element_text(size=12),
        legend.title = element_text(size=14))
citizen_activism_politician_contact_plot
ggsave("citizen_activism_politician_contact_plot.pdf",citizen_activism_politician_contact_plot)
ggsave("citizen_activism_politician_contact_plot.png",citizen_activism_politician_contact_plot)

# FIGURE A6.1 IN SI
# =================

myColorScale <- brewer.pal(5,"Dark2")
citizen_activism_petition_signing_plot <- ggplot(data=subset(df_citizen_activism,activity=="Signed Petition"))+
  geom_pointrange(aes(x=issue_type,
                      y=mean,ymin=cilo,ymax=cihi,group=mean,
                      colour=issue_type,shape=position_on_issue),
                  position=position_dodge(width=0.5),size=0.8)+
  scale_colour_manual(name = "Policy Category",values = myColorScale)+
  scale_shape_manual(name = "Position on Issue",values = c(16,17))+
  xlab("")+
  ylab("Mean Share of Citzens Who Signed a Petition")+
  theme(axis.text.x = element_text(size=12,angle=30,hjust=1),
        axis.text.y = element_text(size=12),
        axis.title.y = element_text(size=12),
        legend.text = element_text(size=12),
        legend.title = element_text(size=14))
citizen_activism_petition_signing_plot
ggsave("citizen_activism_petition_signing_plot.pdf",citizen_activism_petition_signing_plot)
ggsave("citizen_activism_petition_signing_plot.png",citizen_activism_petition_signing_plot)


# DATA TABLE FOR FIGURE A6.1 [SI TABLE A4.4]
# ==========================================

df_citizen_activism

subset(df_citizen_activism,activity=="Signed Petition") %>%
  select(-c(sd,cilo,cihi,activity)) %>%
  as.data.frame() %>%
  stargazer(summary=FALSE,type="html")

